
# reading the data

install.packages("PanelMatch")
library(PanelMatch)
library(readr)
library(haven)

setwd("C:/Users/jttuki/Dropbox/Love thy politician/Results")
rm(list=ls(all=TRUE))

zipcode_sorting2 <- read_dta("C:/Users/jttuki/Dropbox/Love thy politician/data/zipcode_sorting_mainvars.dta")
View(zipcode_sorting2)
summary(zipcode_sorting2)
class(zipcode_sorting2)

#clean data structure.  This read is giving tbl, tbl_df, dataframe types, we want only dataframe
zipcode_sorting <- data.frame(zipcode_sorting2[,1]) 
for (column in c(2:ncol(zipcode_sorting2))){
    print(column)
    zipcode_sorting <- cbind(zipcode_sorting, zipcode_sorting2[,column])
}

class(zipcode_sorting)
class(zipcode_sorting) != "data.frame"
summary(zipcode_sorting)
View(zipcode_sorting)

(class(zipcode_sorting[, year]) != "integer")
class(zipcode_sorting$year)

zipcode_sorting$year<-as.integer(zipcode_sorting$year)
zipcode_sorting$zipcode<-as.integer(zipcode_sorting$zipcode)

class(zipcode_sorting$treataft)
class(zipcode_sorting$log_9_10)
class(zipcode_sorting$opp16)
class(zipcode_sorting$logpop)

estimation_data_small <- zipcode_sorting[zipcode_sorting$opp16<90,]
estimation_data_large <- zipcode_sorting[zipcode_sorting$opp16>89,]
summary(zipcode_sorting)
summary(estimation_data_small)
summary(estimation_data_large)

#Imai et al. 
PM.match.large.high <- PanelMatch(lag=3, time.id= "year", unit.id = "zipcode", treatment = "treataft", refinement.method = "mahalanobis", size.match = 2, covs.formula = ~ opp16 + logpop, data = estimation_data_large, qoi = "att", lead = 0:3, outcome.var = "log_9_10", forbid.treatment.reversal = TRUE, placebo.test = TRUE)
PM.results.large.high <- PanelEstimate(sets = PM.match.large.high, data = estimation_data_large, number.iterations = 100) 
placebo.large.high <- placebo_test(PM.match.large.high,data = estimation_data_large, plot = F, number.iterations = 100)
sdt3lh <- sd(placebo.large.high$bootstrapped.estimates[,1])
sdt2lh <- sd(placebo.large.high$bootstrapped.estimates[,2])
pla_selh <- rbind(sdt3lh,sdt2lh,0)
pla_estlh <- rbind(placebo.large.high$estimates[1],placebo.large.high$estimates[2],0)
pla_rlh<- cbind(pla_estlh,pla_selh)
R.large.high <- cbind(PM.results.large.high$estimates,PM.results.large.high$standard.error)
R.large.high.final<-rbind(pla_rlh,R.large.high)

PM.match.large.mid <- PanelMatch(lag=3, time.id= "year", unit.id = "zipcode", treatment = "treataft", refinement.method = "mahalanobis", size.match = 2, covs.formula = ~ opp16 + logpop, data = estimation_data_large, qoi = "att", lead = 0:3, outcome.var = "log_3_8", forbid.treatment.reversal = TRUE, placebo.test = TRUE)
PM.results.large.mid <- PanelEstimate(sets = PM.match.large.mid, data = estimation_data_large, number.iterations = 100) 
placebo.large.mid <- placebo_test(PM.match.large.mid,data = estimation_data_large, plot = F, number.iterations = 100)
sdt3lm <- sd(placebo.large.mid$bootstrapped.estimates[,1])
sdt2lm <- sd(placebo.large.mid$bootstrapped.estimates[,2])
pla_selm <- rbind(sdt3lm,sdt2lm,0)
pla_estlm <- rbind(placebo.large.mid$estimates[1],placebo.large.mid$estimates[2],0)
pla_rlm<- cbind(pla_estlm,pla_selm)
R.large.mid <- cbind(PM.results.large.mid$estimates,PM.results.large.mid$standard.error)
R.large.mid.final<-rbind(pla_rlm,R.large.mid)

PM.match.large.low <- PanelMatch(lag=3, time.id= "year", unit.id = "zipcode", treatment = "treataft", refinement.method = "mahalanobis", size.match = 2, covs.formula = ~ opp16 + logpop, data = estimation_data_large, qoi = "att", lead = 0:3, outcome.var = "log_1_2", forbid.treatment.reversal = TRUE, placebo.test = TRUE)
PM.results.large.low <- PanelEstimate(sets = PM.match.large.low, data = estimation_data_large, number.iterations = 100) 
placebo.large.low <- placebo_test(PM.match.large.low,data = estimation_data_large, plot = F, number.iterations = 100)
sdt3ll <- sd(placebo.large.low$bootstrapped.estimates[,1])
sdt2ll <- sd(placebo.large.low$bootstrapped.estimates[,2])
pla_sell <- rbind(sdt3ll,sdt2ll,0)
pla_estll <- rbind(placebo.large.low$estimates[1],placebo.large.low$estimates[2],0)
pla_rll<- cbind(pla_estll,pla_sell)
R.large.low <- cbind(PM.results.large.low$estimates,PM.results.large.low$standard.error)
R.large.low.final<-rbind(pla_rll,R.large.low)

PM.match.large.all <- PanelMatch(lag=3, time.id= "year", unit.id = "zipcode", treatment = "treataft", refinement.method = "mahalanobis", size.match = 2, covs.formula = ~ opp16 + logpop, data = estimation_data_large, qoi = "att", lead = 0:3, outcome.var = "logpop", forbid.treatment.reversal = TRUE, placebo.test = TRUE)
PM.results.large.all <- PanelEstimate(sets = PM.match.large.all, data = estimation_data_large, number.iterations = 100) 
placebo.large.all <- placebo_test(PM.match.large.all,data = estimation_data_large, plot = F, number.iterations = 100)
sdt3la <- sd(placebo.large.all$bootstrapped.estimates[,1])
sdt2la <- sd(placebo.large.all$bootstrapped.estimates[,2])
pla_sela <- rbind(sdt3la,sdt2la,0)
pla_estla <- rbind(placebo.large.all$estimates[1],placebo.large.all$estimates[2],0)
pla_rla<- cbind(pla_estla,pla_sela)
R.large.all <- cbind(PM.results.large.all$estimates,PM.results.large.all$standard.error)
R.large.all.final<-rbind(pla_rla,R.large.all)

PM.match.small.high <- PanelMatch(lag=3, time.id= "year", unit.id = "zipcode", treatment = "treataft", refinement.method = "mahalanobis", size.match = 2, covs.formula = ~ opp16 + logpop, data = estimation_data_small, qoi = "att", lead = 0:3, outcome.var = "log_9_10", forbid.treatment.reversal = TRUE, placebo.test = TRUE)
PM.results.small.high <- PanelEstimate(sets = PM.match.small.high, data = estimation_data_small, number.iterations = 100) 
placebo.small.high <- placebo_test(PM.match.small.high,data = estimation_data_small, plot = F, number.iterations = 100)
sdt3sh <- sd(placebo.small.high$bootstrapped.estimates[,1])
sdt2sh <- sd(placebo.small.high$bootstrapped.estimates[,2])
pla_sesh <- rbind(sdt3sh,sdt2sh,0)
pla_estsh <- rbind(placebo.small.high$estimates[1],placebo.small.high$estimates[2],0)
pla_rsh<- cbind(pla_estsh,pla_sesh)
R.small.high <- cbind(PM.results.small.high$estimates,PM.results.small.high$standard.error)
R.small.high.final<-rbind(pla_rsh,R.small.high)

PM.match.small.mid <- PanelMatch(lag=3, time.id= "year", unit.id = "zipcode", treatment = "treataft", refinement.method = "mahalanobis", size.match = 2, covs.formula = ~ opp16 + logpop, data = estimation_data_small, qoi = "att", lead = 0:3, outcome.var = "log_3_8", forbid.treatment.reversal = TRUE, placebo.test = TRUE)
PM.results.small.mid <- PanelEstimate(sets = PM.match.small.mid, data = estimation_data_small, number.iterations = 100) 
placebo.small.mid <- placebo_test(PM.match.small.mid,data = estimation_data_small, plot = F, number.iterations = 100)
sdt3sm <- sd(placebo.small.mid$bootstrapped.estimates[,1])
sdt2sm <- sd(placebo.small.mid$bootstrapped.estimates[,2])
pla_sesm <- rbind(sdt3sm,sdt2sm,0)
pla_estsm <- rbind(placebo.small.mid$estimates[1],placebo.small.mid$estimates[2],0)
pla_rsm<- cbind(pla_estsm,pla_sesm)
R.small.mid <- cbind(PM.results.small.mid$estimates,PM.results.small.mid$standard.error)
R.small.mid.final<-rbind(pla_rsm,R.small.mid)

PM.match.small.low <- PanelMatch(lag=3, time.id= "year", unit.id = "zipcode", treatment = "treataft", refinement.method = "mahalanobis", size.match = 2, covs.formula = ~ opp16 + logpop, data = estimation_data_small, qoi = "att", lead = 0:3, outcome.var = "log_1_2", forbid.treatment.reversal = TRUE, placebo.test = TRUE)
PM.results.small.low <- PanelEstimate(sets = PM.match.small.low, data = estimation_data_small, number.iterations = 100) 
placebo.small.low <- placebo_test(PM.match.small.low,data = estimation_data_small, plot = F, number.iterations = 100)
sdt3sl <- sd(placebo.small.low$bootstrapped.estimates[,1])
sdt2sl <- sd(placebo.small.low$bootstrapped.estimates[,2])
pla_sesl <- rbind(sdt3sl,sdt2sl,0)
pla_estsl <- rbind(placebo.small.low$estimates[1],placebo.small.low$estimates[2],0)
pla_rsl<- cbind(pla_estsl,pla_sesl)
R.small.low <- cbind(PM.results.small.low$estimates,PM.results.small.low$standard.error)
R.small.low.final<-rbind(pla_rsl,R.small.low)

PM.match.small.all <- PanelMatch(lag=3, time.id= "year", unit.id = "zipcode", treatment = "treataft", refinement.method = "mahalanobis", size.match = 2, covs.formula = ~ opp16 + logpop, data = estimation_data_small, qoi = "att", lead = 0:3, outcome.var = "logpop", forbid.treatment.reversal = TRUE, placebo.test = TRUE)
PM.results.small.all <- PanelEstimate(sets = PM.match.small.all, data = estimation_data_small, number.iterations = 100) 
placebo.small.all <- placebo_test(PM.match.small.all,data = estimation_data_small, plot = F, number.iterations = 100)
sdt3sa <- sd(placebo.small.all$bootstrapped.estimates[,1])
sdt2sa <- sd(placebo.small.all$bootstrapped.estimates[,2])
pla_sesa <- rbind(sdt3sa,sdt2sa,0)
pla_estsa <- rbind(placebo.small.all$estimates[1],placebo.small.all$estimates[2],0)
pla_rsa<- cbind(pla_estsa,pla_sesa)
R.small.all <- cbind(PM.results.small.all$estimates,PM.results.small.all$standard.error)
R.small.all.final<-rbind(pla_rsa,R.small.all)

x<-c(-3,-2,-1,0,1,2,3)
time<-c("t-3","t-2","t-1","t+0","t+1","t+2","t+3")

ImaiData <- cbind(x,time,R.large.high.final,R.large.mid.final,R.large.low.final,R.large.all.final,R.small.high.final,R.small.mid.final,R.small.low.final,R.small.all.final)
ImaiData2 <- data.frame(ImaiData[,1]) 
for (column in c(2:ncol(ImaiData))){
  print(column)
  ImaiData2 <- cbind(ImaiData2, ImaiData[,column])
}
colnames(ImaiData2) <- c("x","time","est_l_hi","sd_l_hi","est_l_mi","sd_l_mi","est_l_lo","sd_l_lo","est_l_all","sd_l_all","est_s_hi","sd_s_hi","est_s_mi","sd_s_mi","est_s_lo","sd_s_lo","est_s_all","sd_s_all") 
write_dta(ImaiData2,"C:/Users/jttuki/Dropbox/Love thy politician/dofiles/Journal/JOP_replication/DID_est_Imai_et_al_replication.dta")

